home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Language/OS - Multiplatform Resource Library
/
LANGUAGE OS.iso
/
cpp_libs
/
rwvector.lha
/
RWVector2.1
/
src
/
testdfft.cc
< prev
next >
Wrap
C/C++ Source or Header
|
1989-08-18
|
4KB
|
140 lines
/*
* Test real FFT server
*
* Copyright (C) 1988, 1989.
*
* Dr. Thomas Keffer
* Rogue Wave Associates
* P.O. Box 85341
* Seattle WA 98145-1341
*
* Permission to use, copy, modify, and distribute this
* software and its documentation for any purpose and
* without fee is hereby granted, provided that the
* above copyright notice appear in all copies and that
* both that copyright notice and this permission notice
* appear in supporting documentation.
*
* This software is provided "as is" without any
* expressed or implied warranty.
*
*
* @(#)testdfft.cc 2.1 8/18/89
*/
#include "rw/DoubleFFT.h"
overload printVec;
static void
printVec(ostream& s, const DComplexVec& x);
static void
printVec(ostream& s, const DoubleVec& x);
main()
{
// npts must be even:
const unsigned npts = 12;
cout <<"Testing DoubleFFTServer (Double Precision FFT Server)\n\n";
DoubleFFTServer server;
DoubleVec a(cos(DoubleVec(npts, 0, 2.0*M_PI/npts))); // One period
DoubleVec a2(sin(DoubleVec(npts, 0, 4.0*M_PI/npts))); // Two periods
DoubleVec asum = a + a2;
cout <<"**************************************\n";
cout <<"a:\n";
printVec(cout,a);
cout <<"a2:\n";
printVec(cout,a2);
cout <<"asum:\n";
printVec(cout, asum);
cout <<"**************************************\n";
cout <<"Checking transforms of pure signals.\n";
cout <<"\nTransform of a:\n";
printVec(cout,server.fourier(a)/DComplex(npts));
cout <<"\nTransform of a2:\n";
printVec(cout, server.fourier(a2)/DComplex(npts));
cout <<"\nTransform of asum:\n";
// asum_fourier will be a complex conjugate even sequence
DComplexVec asum_fourier = server.fourier(asum)/DComplex(npts);
printVec(cout, asum_fourier);
cout <<"**************************************\n";
cout <<"Checking Parseval's theorem.\n";
cout <<"\nOriginal variance: "<<variance(asum)<<NL;
// expand to full complex series:
double var = spectralVariance(expandConjugateEven(asum_fourier));
cout <<"Spectral variance: "<<var<<NL;
cout <<"**************************************\n";
DoubleVec aramp(npts,0,1);
cout <<"Checking transform of linear ramp.\n";
cout <<"\nOriginal series:\n";
printVec(cout, aramp);
DComplexVec aramp_fourier = server.fourier(aramp)/DComplex(npts);
cout <<"Its transform:\n";
printVec(cout, aramp_fourier);
cout <<NL;
cout <<"Checking Parseval's theorem for ramp.\n";
cout <<"Original variance: "<<variance(aramp)<<NL;
var = spectralVariance(expandConjugateEven(aramp_fourier));
cout <<"Final variance: "<<var<<NL;
cout <<"\nBack transform:\n";
printVec(cout, server.ifourier(aramp_fourier) );
cout <<"**************************************\n";
cout <<"Checking Nyquist.\n";
DoubleVec Nyquist(npts,1.0);
Nyquist.slice(1,npts/2,2) = -1.0;
cout <<"Original sequence:\n";
printVec(cout, Nyquist);
cout <<"Its transform:\n";
DComplexVec Nyquist_fourier = server.fourier(Nyquist)/DComplex(npts);
printVec(cout, Nyquist_fourier);
cout <<"Back transform:\n";
printVec(cout, server.ifourier(Nyquist_fourier));
cout <<"**************************************\n";
exit(0);
}
// Use special functions to print vectors, to round the output
// to a fixed point number, so that differences at the machine
// precision won't be seen.
static void
printVec(ostream& s, const DComplexVec& x)
{
for(int i = 0; i<x.length(); i++){
double rl = real(x(i));
double img = imag(x(i));
if(!(i%4) && i) s<<"\n";
s << "(" << form("%8.5f", rl) << ", ";
s << form("%8.5f", img) << ")" ;
}
s<<"\n";
}
static void
printVec(ostream& s, const DoubleVec& x)
{
for(int i = 0; i<x.length(); i++){
if(!(i%5) && i) s<<"\n";
s << form("%8.5f",x(i)) << " ";
}
s<<"\n";
}